#!/usr/bin/env python3
# -*- coding: utf-8 -*-
# Author: leuelcey (wyd)
# E-mail: leuelcey@163.com
# Creation date: 2021.2
# Modified date: 2021.5

from math import sin, cos, tan, asin, acos, atan, atan2, floor, pi as PI
from functools import lru_cache
from basicdata import BasicData
tbd = BasicData()
Radian = tbd.get_constant(section='astronomy', varname='radian_arcsec')
J2000 = tbd.get_constant(section='astronomy', varname='J2000')
Pi2 = tbd.get_constant(section='astronomy', varname='pi_2times')

# 基准儒略日
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def base_Juday(bY, bM):
    return Gregory2Juday(bY, bM, 1 + ((12 + 0.1 / 3600) / 24))

# 二次曲线外推
# dt_ext, conic_extrapolation
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def conic_extrapolate(y, jsd):
    return -20 + jsd * ((y - 1820) / 100) * ((y - 1820) / 100)

# 计算世界时与原子时之差(TD-UT，单位:日)，y为西元年
# dt_calc, ut_at_Dvalue
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def UTAT_Dvalue(y):
    dt_at = tbd.get_constant(section='deltaT')
    y0 = dt_at[len(dt_at) - 2]
    if y >= y0:
        if y > y0 + 100:
            return conic_extrapolate(y, 31)
        dv = conic_extrapolate(y0, 31) - dt_at[len(dt_at) - 1]
        return conic_extrapolate(y, 31) - dv * (y0 + 100 - y) / 100

    for i in range(0, len(dt_at), 5):
        if y < dt_at[i + 5]:
            break
    t1 = (y - dt_at[i]) / (dt_at[i + 5] - dt_at[i]) * 10
    return dt_at[i + 1] + dt_at[i + 2] * t1 + dt_at[i + 3] * t1 * t1 + dt_at[i + 4] * t1 * t1 * t1

# delta T计算，t为儒略日(J2000起算)
# dt_T
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def deltaT(t):
    return UTAT_Dvalue(t / 365.2425 + 2000.0) / 86400.0

# 除冥王星外的八大行星的星历计算，xt为星体序号，zn为坐标号，t为儒略世纪数，n为计算项数
# XL0_calc
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
@lru_cache
def planet_ephemeris(xt, zn, t, n):
    t = t / 10
    v, tn = 0, 1

    F = tbd.get_constant(section='Ephemeris', varname='Planets', subpart=xt)
    pn = zn * 6 + 1
    N0 = F[pn + 1] - F[pn]
    for i in range(6):
        n1 = F[pn + i]
        n2 = F[pn + 1 + i]
        n0 = n2 - n1
        if not n0:
            continue
        if n < 0:
            N = n2
        else:
            N = floor(3 * n * n0 / N0 + 0.5) + n1
            if i:
                N = N + 3
            if N > n2:
                N = n2
        c = 0
        for j in range(n1, N, 3):
            c = c + (F[j] * cos(F[j + 1] + t * F[j + 2]))
        v = v + c * tn
        tn = tn * t

    v = v / F[0]
    if xt == 0:
        if zn == 0:
            v = v + ((-0.0728 - 2.7702 * t - 1.1019 * t * t - 0.0996 * t * t * t) / Radian)
        if zn == 1:
            v = v + ((+0.0000 + 0.0004 * t + 0.0004 * t * t - 0.0026 * t * t * t) / Radian)
        if zn == 2:
            v = v + ((-0.0020 + 0.0044 * t + 0.0213 * t * t - 0.0250 * t * t * t) / 1000000)
    else:
        if zn == 0:
            v = v + (-3 * t / Radian)
        if zn == 2:
            dv = tbd.get_constant(section='Ephemeris', varname='Correct')[(xt - 1) * 3 + zn]
            v = v + dv / 1000000
        else:
            dv = tbd.get_constant(section='Ephemeris', varname='Correct')[(xt - 1) * 3 + zn]
            v = v + dv / Radian
    return v

# 地球经度计算，返回Date分点黄经，传入世纪数和取项数
# E_Lon
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def earth_lng(t, n):
    return planet_ephemeris(0, 0, t, n)

# 黄经章动
# nutationLon2
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
@lru_cache
def nutation_in_lng(t):
    t2 = t * t
    dL = 0
    B = tbd.get_constant(section='Nutation', varname='Medium')
    for i in range(0, len(B), 5):
        if i == 0:
            a = -1.742 * t
        else:
            a = 0
        dL = dL + ((B[i + 3] + a) * sin(B[i] + B[i + 1] * t + B[i + 2] * t2))
    return dL / 100 / Radian

# 太阳光行差，t为世纪数
# gxc_sunLon
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def sun_aberration(t):
    v = -0.043126 + 628.301955 * t - 0.000002732 * t * t
    e = 0.016708634 - 0.000042037 * t - 0.0000001267 * t * t
    return (-20.49552 * (1 + e * cos(v))) / Radian

# 太阳视黄经
# S_aLon, sun_apparent_lng
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def sun_alng(t, n):
    return earth_lng(t, n) + nutation_in_lng(t) + sun_aberration(t) + PI

# 对超过0-2PI的角度转为0-2PI
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def rad2mrad(v):
    v = v % (2 * PI)
    if v < 0:
        return v + 2 * PI
    return v

# 球面坐标旋转，黄道赤道坐标变换(赤到黄时E取负)
# llrConv
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def coordinate_rotation(JW, E):
    r = [0.0, 0.0, 0.0]
    J = JW[0]
    W = JW[1]
    r[0] = atan2(sin(J) * cos(E) - tan(W) * sin(E), cos(J))
    r[1] = asin(cos(E) * sin(W) + sin(E) * cos(W) * sin(J))
    r[2] = JW[2]
    r[0] = rad2mrad(r[0])
    return r

# 对超过-PI到PI的角度转为-PI到PI
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def rad2rrad(v):
    v = v % (2 * PI)
    if v <= -PI:
        return v + 2 * PI
    if v > PI:
        return v - 2 * PI
    return v

# 平太阳真太阳时差计算(低精度)，误差约在1秒以内，t为力学时儒略世纪数
# pty_zty2
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def equation_of_time(t):
    L = (1753470142 + 628331965331.8 * t + 5296.74 * t * t) / 1000000000 + PI
    z = [0.0, 0.0, 0.0]
    E = (84381.4088 - 46.836051 * t) / Radian
    z[0] = planet_ephemeris(0, 0, t, 5) + PI
    z = coordinate_rotation(z, E)
    L = rad2rrad(L - z[0])
    return L / Pi2

# 儒略日转为时间串(hh:mm:ss)
# timeStr
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def Juday2timestring(jd):
    jd = jd + 0.5
    jd = (jd - floor(jd))
    s = floor(jd * 86400 + 0.5)
    h = floor(s / 3600)
    s = s - (h * 3600)
    m = floor(s / 60)
    s = s - (m * 60)
    sh = "0" + str(h)
    sm = "0" + str(m)
    ss = "0" + str(s)
    return sh[-2:] + ':' + sm[-2:] + ':' + ss[-2:]

# 专用，时间串(hh:mm:ss)转换为小时值
# timeStr2hour
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def timestring2hour(s):
    s = s.strip().split(':')
    hh, mm, ss = 0, 0, 0
    if len(s) == 1:
        hh = int(s[0])
        mm = int(s[0][2:4])
        ss = int(s[0][4:2])
    if len(s) == 2:
        hh = int(s[0])
        mm = int(s[1])
        ss = 0
    if len(s) == 3:
        hh = int(s[0])
        mm = int(s[1])
        ss = int(s[2])
    return hh + mm / 60 + ss / 3600

# 本地平太阳时计算。secs为本日过去的秒数,J为本地经度
# bdptys
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def mean_solar_time(secs, J, Gs=''):
    bzjd0 = 120  # 标准时间发出地经度
    ptys = floor(secs + 4 * (J - bzjd0) * 60)  # 平太阳时,每经度时差4分钟
    ptysds = ptys / 3600 / 24  # 平太阳时化为日数

    if ptys < 0:
        ptys = 24 * 3600 + ptys  # 负值时转为前一日
    ptysh = floor(ptys / 3600)  # 平小时
    ptysm = floor((ptys - ptysh * 3600) / 60)  # 平分钟
    ptyss = floor(ptys - ptysh * 3600 - ptysm * 60)  # 平秒
    sh = '0' + str(ptysh) if (ptysh < 10) else str(ptysh)
    sm = '0' + str(ptysm) if (ptysm < 10) else str(ptysm)
    ss = '0' + str(ptyss) if (ptyss < 10) else str(ptyss)

    if Gs == 'd':
        return ptysds  # 返回平太阳日数
    else:
        return sh + ':' + sm + ':' + ss  # 返回平太阳时串

# 专用，将输入的年份进行通用化处理，返回天文纪年
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def year2Ayear(c):
    y = 0
    c = c.strip()
    q = c[0]
    if q in ('B', 'b', '-'):
        y = 1 - int(c[1:])
        if y > 0:
            #通用纪法的公元前纪法从B.C.1年开始。并且没有公元0年
            return -10000
    else:
        y = int(c)

    if y < -4712:
        # 超过B.C. 4713不准
        y = -4712
    if y > 9999:
        # 超过9999年的农历计算很不准
        y = 9999
    return y

# 西历时间转换为儒略日数，参数为西元年月日，返回对应的儒略日数
# JD
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
@lru_cache
def Gregory2Juday(y, m, d):
    n = 0
    G = 0
    if (y * 372 + m * 31 + floor(d)) >= 588829:
        G = 1
    if m <= 2:
        m = m + 12
        y = y - 1
    if G > 0:
        n = floor(y / 100)
        n = 2 - n + floor(n / 4)
    return floor(365.25 * (y + 4716)) + floor(30.6001 * (m + 1)) + d + n - 1524.5

# 儒略日数转换为西历时间，jd 为儒略日数，返回西元年月日时分秒
# DD
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
@lru_cache
def Juday2Gregory(jd):
    D = floor(jd + 0.5)
    F = jd + 0.5 - D
    if D >= 2299161:
        c = floor((D - 1867216.25) / 36524.25)
        D = D + (1 + c - floor(c / 4))

    D = D + 1524
    xY = floor((D - 122.1) / 365.25)
    D = D - floor(365.25 * xY)
    xM = floor(D / 30.601)
    D = D - floor(30.601 * xM)
    xD = D
    if xM > 13:
        xM = xM - 13
        xY = xY - 4715
    else:
        xM = xM - 1
        xY = xY - 4716

    F = F * 24
    xh = floor(F)
    F = F - xh
    F = F * 60
    xm = floor(F)
    F = F - xm
    F = F * 60
    xs = F
    return xY, xM, xD, xh, xm, round(xs, 3)

# 专用，经纬度解码，返回经度和维度值
# JWdecode
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def lng_lat_decode(v):
    jw = [0, 0, 0, 0]
    for i in range(4):
        jw[i] = ord(v[i])
        if jw[i] > 96:
            jw[i] = jw[i] - (97 - 36)
        elif jw[i] > 64:
            jw[i] = jw[i] - (65 - 10)
        else:
            jw[i] = jw[i] - 48

    lsJ = (jw[2] + jw[3] / 60 + 73) / 180 * PI
    lsW = (jw[0] + jw[1] / 60) / 180 * PI
    return round(lsJ / PI * 180, 6), round(lsW / PI * 180, 6)

# 根据经度获取时区值，J 为经度值，返回时区值
# TZcalc
# by leuelcey (wyd)
def get_timezone(J):
    shang, yushu = divmod(J, 15)  # 取得经度与15的商和余数
    if yushu < 7.5:  # 若余数小于7.5，则时区值即为商数
        return shang
    if yushu > 7.5:  # 若余数大于7.5，则时区值为商数加1
        return shang + 1
    return 0

# 地球速度，t为世纪数，误差小于万分3
# E_v
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def earth_velocity(t):
    f = 628.307585 * t
    return 628.332 + 21 * sin(1.527 + f) + 0.44 * sin(1.48 + f * 2) + 0.129 * sin(5.82 + f) * t + 0.00055 * sin(4.21 + f) * t * t

# 已知太阳视黄经反求时间
# S_aLon_t, apparent_lng2time
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def sun_alng2time(W):
    v = 628.3319653318
    t = (W - 1.75347 - PI) / v
    v = earth_velocity(t)
    t = t + (W - sun_alng(t, 10)) / v
    v = earth_velocity(t);
    t = t + (W - sun_alng(t, -1)) / v
    return t

# 月亮星历计算
# XL1_calc
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
@lru_cache
def moon_ephemeris(zn, t, n):
    ob = tbd.get_constant(section='Ephemeris', varname='Moon', subpart=zn)
    v = 0
    tn = 1
    t2 = t * t
    t3 = t2 * t
    t4 = t3 * t
    t5 = t4 * t
    tx = t - 10
    if zn == 0:
        v = v + ((3.81034409 + 8399.684730072 * t - 3.319e-05 * t2 + 3.11e-08 * t3 - 2.033e-10 * t4) * Radian)
        v = v + (5028.792262 * t + 1.1124406 * t2 + 0.00007699 * t3 - 0.000023479 * t4 - 0.0000000178 * t5)
        if tx > 0:
            v = v + (-0.866 + 1.43 * tx + 0.054 * tx * tx)

    t2 = t2 / 1e4
    t3 = t3 / 1e8
    t4 = t4 / 1e8
    n = n * 6
    if n < 0:
        n = len(ob[0])
    for i in range(len(ob)):
        F = ob[i]
        N = floor(n * len(F) / len(ob[0]) + 0.5)
        if i != 0:
            N = N + 6
        if N >= len(F):
            N = len(F)
        c = 0
        for j in range(0, N, 6):
            c = c + (F[j] * cos(F[j + 1] + t * F[j + 2] + t2 * F[j + 3] + t3 * F[j + 4] + t4 * F[j + 5]))
        v = v + (c * tn)
        tn = tn * t

    if zn != 2:
        v = v / Radian
    return v

# 月球经度计算，返回Date分点黄经，t为世纪数，n为项数比例
# M_Lon
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def moon_Lng(t, n):
    return moon_ephemeris(0, t, n)

# 月球经度光行差，误差0.07"
# gxc_moonLon
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def moon_aberration(t):
    return -3.4E-6

# 月日视黄经的差值
# MS_aLon, moonsun_lng_Dvalue
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def ms_alng_Dvalue(t, Mn, Sn):
    return moon_Lng(t, Mn) + moon_aberration(t) - (earth_lng(t, Sn) + sun_aberration(t) + PI)

# 月球速度计算，t为世经数
# M_v
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def moon_velocity(t):
    v = 8399.71 - 914 * sin(0.7848 + 8328.691425 * t + 0.0001523 * t * t)
    v = v - (179 * sin(2.543 + 15542.7543 * t) + 160 * sin(0.1874 + 7214.0629 * t) + 62 * sin(3.14 + 16657.3828 * t) + 34 * sin(4.827 + 16866.9323 * t) + 22 * sin(4.9 + 23871.4457 * t) + 12 * sin(2.59 + 14914.4523 * t) + 7 * sin(0.23 + 6585.7609 * t) + 5 * sin(0.9 + 25195.624 * t) + 5 * sin(2.32 - 7700.3895 * t) + 5 * sin(3.88 + 8956.9934 * t) + 5 * sin(0.49 + 7771.3771 * t))
    return v

# 已知月日视黄经差求时间
# MS_aLon_t, mslng_Dvalue2time
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
def ms_alng2time(W):
    t = 0
    v = 7771.37714500204
    t = (W + 1.08472) / v
    t = t + ((W - ms_alng_Dvalue(t, 3, 3)) / v)
    v = moon_velocity(t) - earth_velocity(t)
    t = t + ((W - ms_alng_Dvalue(t, 20, 10)) / v)
    t = t + ((W - ms_alng_Dvalue(t, -1, 60)) / v)
    return t

# 已知太阳视黄经反求时间,高速低精度,最大误差不超过600秒
# S_aLon_t2, apparent_lng2time2
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
@lru_cache
def sun_alng2time2(W):
    v = 628.3319653318
    t = (W - 1.75347 - PI) / v
    t = t - ((0.000005297 * t * t + 0.0334166 * cos(4.669257 + 628.307585 * t) + 0.0002061 * cos(2.67823 + 628.307585 * t) * t) / v)
    t = t + ((W - earth_lng(t, 8) - PI + (20.5 + 17.2 * sin(2.1824 - 33.75705 * t)) / Radian) / v)
    return t

# 已知月日视黄经差求时间,高速低精度,误差不超过600秒(只验算了几千年)
# MS_aLon_t2, mslng_Dvalue2time2
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
@lru_cache
def ms_alng2time2(W):
    v = 7771.37714500204
    t = (W + 1.08472) / v
    t2 = t * t
    t = t - ((-0.00003309 * t2 + 0.10976 * cos(0.784758 + 8328.6914246 * t + 0.000152292 * t2) + 0.02224 * cos(0.18740 + 7214.0628654 * t - 0.00021848 * t2) - 0.03342 * cos(4.669257 + 628.307585 * t)) / v)
    L = moon_Lng(t, 20) - (4.8950632 + 628.3319653318 * t + 0.000005297 * t * t + 0.0334166 * cos(4.669257 + 628.307585 * t) + 0.0002061 * cos(2.67823 + 628.307585 * t) * t + 0.000349 * cos(4.6261 + 1256.61517 * t) - 20.5 / Radian)
    v = 7771.38 - 914 * sin(0.7848 + 8328.691425 * t + 0.0001523 * t * t) - 179 * sin(2.543 + 15542.7543 * t) - 160 * sin(0.1874 + 7214.0629 * t)
    t = t + ((W - L) / v)
    return t

# 实气实朔计算，用于计算节气、朔日等
# 算法及源代码：许剑伟（莆田十中）
# Python 代码：leuelcey (wyd)
class Climoon(object):
    """
    气朔计算类。
    """
    __slots__ = ('SB', 'QB', 'ZQ', 'HS', 'leap', 'ym', 'dx', 'Yn', 'pe')

    def __init__(self):
        super(Climoon, self).__init__()
        self.SB = ''
        self.QB = ''
        self.leap = 0
        self.ym = []
        self.ZQ = []
        self.HS = []
        self.dx = []
        self.Yn = []
        self.pe = {}
        self._amend_decode()

    @lru_cache
    def so_low(self, W):
        v = 7771.37714500204
        t = (W + 1.08472) / v
        t = t - ((-0.0000331 * t * t + 0.10976 * cos(0.785 + 8328.6914 * t) + 0.02224 * cos(0.187 + 7214.0629 * t) - 0.03342 * cos(4.669 + 628.3076 * t)) / v + (32 * (t + 1.8) * (t + 1.8) - 20) / 86400 / 36525)
        return t * 36525 + 8 / 24

    @lru_cache
    def qi_low(self, W):
        v = 628.3319653318
        t = (W - 4.895062166) / v
        t = t - ((53 * t * t + 334116 * cos(4.67 + 628.307585 * t) + 2061 * cos(2.678 + 628.3076 * t) * t) / v / 10000000)
        L = 48950621.66 + 6283319653.318 * t + 53 * t * t + 334166 * cos(4.669257 + 628.307585 * t) + 3489 * cos(4.6261 + 1256.61517 * t) + 2060.6 * cos(2.67823 + 628.307585 * t) * t - 994 - 834 * sin(2.1824 - 33.75705 * t)
        t = t - ((L / 10000000 - W) / 628.332 + (32 * (t + 1.8) * (t + 1.8) - 20) / 86400 / 36525)
        return t * 36525 + 8 / 24

    @lru_cache
    def qi_high(self, W):
        t = sun_alng2time2(W) * 36525
        t = t - deltaT(t) + 8 / 24
        v = ((t + 0.5) % 1) * 86400
        if v < 1200 or v > 86400 - 1200:
            t = sun_alng2time2(W) * 36525 - deltaT(t) + 8 / 24
        return t

    @lru_cache
    def so_high(self, W):
        t = ms_alng2time2(W) * 36525
        t = t - deltaT(t) + 8 / 24
        v = ((t + 0.5) % 1) * 86400
        if v < 1800 or v > 86400 - 1800:
            t = ms_alng2time(W) * 36525 - deltaT(t) + 8 / 24
        return t

    @lru_cache
    def _amend_decode(self):
        def _practical(s):
            o = "0000000000"
            s = s.replace('J', '00')
            s = s.replace('I', '000')
            s = s.replace('H', '0000')
            s = s.replace('G', '00000')
            s = s.replace('t', '02')
            s = s.replace('s', '002')
            s = s.replace('r', '0002')
            s = s.replace('q', '00002')
            s = s.replace('p', '000002')
            s = s.replace('o', '0000002')
            s = s.replace('n', '00000002')
            s = s.replace('m', '000000002')
            s = s.replace('l', '0000000002')
            s = s.replace('k', '01')
            s = s.replace('j', '0101')
            s = s.replace('i', '001')
            s = s.replace('h', '001001')
            s = s.replace('g', '0001')
            s = s.replace('f', '00001')
            s = s.replace('e', '000001')
            s = s.replace('d', '0000001')
            s = s.replace('c', '00000001')
            s = s.replace('b', '000000001')
            s = s.replace('a', '0000000001')
            s = s.replace('A', o * 6)
            s = s.replace('B', o * 5)
            s = s.replace('C', o * 4)
            s = s.replace('D', o * 3)
            s = s.replace('E', o * 2)
            s = s.replace('F', o)
            return s

        suoS = tbd.get_constant(section='Annles', varname='NewmoonAmend')
        qiS = tbd.get_constant(section='Annles', varname='ClimateAmend')
        self.SB = _practical(suoS)
        self.QB = _practical(qiS)

    @lru_cache
    def qs_Juday(self, jd, qs):
        jd = jd + 2451545
        B = tbd.get_constant(section='Annles', varname='NewmoonLfit')
        pc = 14
        if qs == '气':
            B = tbd.get_constant(section='Annles', varname='ClimateLfit')
            pc = 7

        f1 = int(B[0]) - pc
        f2 = int(B[-1]) - pc
        f3 = 2436935
        if jd < f1 or jd >= f3:
            if qs == '气':
                return floor(self.qi_high(floor((jd + pc - 2451259) / 365.2422 * 24) * PI / 12) + 0.5)
            return floor(self.so_high(floor((jd + pc - 2451551) / 29.5306) * PI * 2) + 0.5)

        if f1 <= jd < f2:
            for i in range(0, len(B), 2):
                if jd + pc < B[i + 2]:
                    break
            D = B[i] + B[i + 1] * floor((jd + pc - B[i]) / B[i + 1])
            D = floor(D + 0.5)
            if D == 1683460:
                D = D + 1
            return D - 2451545

        if f2 <= jd < f3:
            if qs == '气':
                D = floor(self.qi_low(floor((jd + pc - 2451259) / 365.2422 * 24) * PI / 12) + 0.5)
                n = self.QB[floor((jd - f2) / 365.2422 * 24)]
            else:
                D = floor(self.so_low(floor((jd + pc - 2451551) / 29.5306) * PI * 2) + 0.5)
                n = self.SB[floor((jd - f2) / 29.5306)]

            if n == "1":
                return D + 1
            if n == "2":
                return D - 1
            return D

    @lru_cache
    def yue_Sort(self, jd):
        W = floor((jd - 355 + 183) / 365.2422) * 365.2422 + 355
        if self.qs_Juday(W, '气') > jd:
            W = W - 365.2422

        if not self.ZQ:
            self.ZQ = [0 for x in range(25)]
        for i in range(25):
            self.ZQ[i] = self.qs_Juday(W + 15.2184 * i, '气')

        self.pe['pe1'] = self.qs_Juday(W - 15.2, '气')
        self.pe['pe2'] = self.qs_Juday(W - 30.4, '气')
        w = self.qs_Juday(self.ZQ[0], '朔')
        if w > self.ZQ[0]:
            w = w - 29.53

        if not self.HS:
            self.HS = [0 for x in range(15)]
        for i in range(15):
            self.HS[i] = self.qs_Juday(w + 29.5306 * i, '朔')

        self.leap = 0
        if not self.dx:
            self.dx = [0 for x in range(14)]
        nym = [0 for x in range(14)]
        for i in range(14):
            self.dx[i] = self.HS[i + 1] - self.HS[i]
            nym[i] = i

        ymc = tbd.get_constant(section='calendar', varname='yming')
        self.ym = ['' for x in range(14)]
        YY = floor((self.ZQ[0] + 10 + 180) / 365.2422) + 2000

        if -721 <= YY <= -104:
            ns = [0, 0, 0, '', '', '', 0, 0, 0]
            yy = 0
            for i in range(3):
                yy = YY + i - 1
                if yy >= -721:
                    ns[i] = self.qs_Juday(1457698 - J2000 + floor(0.342 + (yy + 721) * 12.368422) * 29.5306, '朔')
                    ns[i + 3] = '十三'
                    ns[i + 6] = 2
                if yy >= -479:
                    ns[i] = self.qs_Juday(1546083 - J2000 + floor(0.500 + (yy + 479) * 12.368422) * 29.5306, '朔')
                    ns[i + 3] = '十三'
                    ns[i + 6] = 2
                if yy >= -220:
                    ns[i] = self.qs_Juday(1640641 - J2000 + floor(0.866 + (yy + 220) * 12.369000) * 29.5306, '朔')
                    ns[i + 3] = '后九'
                    ns[i + 6] = 11

            nn = 0
            for i in range(14):
                for nn in range(2, -1, -1):
                    if self.HS[i] >= ns[nn]:
                        break
                f1 = floor((self.HS[i] - ns[nn] + 15) / 29.5306)
                if f1 < 12:
                    self.ym[i] = ymc[(f1 + ns[nn + 6]) % 12]
                else:
                    self.ym[i] = ns[nn + 3]
            return

        if self.HS[13] <= self.ZQ[24]:
            i = 1
            while self.HS[i + 1] > self.ZQ[2 * i] and i < 13:
                i = i + 1
            self.leap = i
            for j in range(i, 14):
                nym[j] = nym[j] - 1

        for i in range(14):
            Dm = self.HS[i] + J2000
            v2 = nym[i]
            mc = ymc[v2 % 12]
            if (1724360 <= Dm <= 1729794) or (1807724 <= Dm <= 1808699):
                mc = ymc[(v2 + 1) % 12]
            if 1999349 <= Dm <= 1999467:
                mc = ymc[(v2 + 2) % 12]
            if 1973067 <= Dm <= 1977052:
                if v2 % 12 == 0:
                    mc = "正"
                if v2 == 2:
                    mc = '一'

            if Dm == 1729794 or Dm == 1808699:
                mc = '拾贰'
            self.ym[i] = mc

# by leuelcey (wyd)
class sVar(object):
    """
    用于模拟静态变量。
    """
    __datas = None

    @classmethod
    def Dump(cls, flag=None, data=None):
        if flag is not None and data is not None:
            if not cls.__datas:
                cls.__datas = {}
            if flag in cls.__datas:
                if cls.__datas[flag] != data:
                    cls.__datas[flag] = data
            else:
                cls.__datas[flag] = data

    @classmethod
    def Load(cls, flag=None):
        if flag is not None and cls.__datas is not None:
            if flag in cls.__datas:
                return cls.__datas[flag]
            return None
        return None

    @classmethod
    def List(cls):
        if cls.__datas:
            for _, flag in enumerate(cls.__datas):
                print(flag)

    @classmethod
    def Init(cls):
        if cls.__datas is not None:
            cls.__datas = None

# 完成类似于 PHP 的 array_keys 的功能
# by leuelcey (wyd)
def get_keys(Array=None, Value=None, Default=None):
    if not Array or not isinstance(Array, (list, tuple)):
        return []

    res = []
    for k, v in enumerate(Array):
        if v == Value:
            res = res + [k if not Default else Default]
    return res

